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Abstract 



We calculate the leading and next-to-leading logarithmic resummed distri- 
bution for the jet broadening in deep inelastic scattering, as well as the 
power correction for both the distribution and mean value. A truncation 
of the answer at NLL accuracy, as is standard, leads to unphysical diver- 
gences. We discuss their origin and show how the problem can be resolved. 
We then examine DIS-specific procedures for matching to fixed-order calcu- 
lations and compare our results to data. One of the tools developed for the 
comparison is an NLO parton distribution evolution code. When compared 
to PDF sets from MRST and CTEQ it reveals limited discrepancies in both. 
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1 Introduction 



In e~^e~ collisions there have been extensive studies of event shape distributions involv- 
ing comparisons to calculations which resum logarithms at the edge of phase space . 
Much has been learnt from these studies, for example precise determinations of as and 
strong tests of recent novel approaches to hadronisation (see for example f^), and even 
explicit measurements of the colour factors of QCD [0. 

In the past few years the HI and ZEUS experiments at HERA have embarked on 
analogous studies of DIS current-hemisphere event shapes [PHTOH- Compared to most 



e~^e~ results, a feature of the DIS measurements is that a wide range of Q values is 
probed by the same experiment. This is useful when one wishes to isolate effects with 
a specific dependence on Q, such as the running of or hadronisation corrections. 
Additionally DIS event shapes depend on radiation from the incoming proton, allowing 
one to study non-trivial questions related to the process-independence of hadronisation 
corrections and perhaps even issues such as intrinsic transverse momentum. 

With these motivations in mind we recently initiated a project to resum a range of 



event shapes in DIS |11-13]. This paper deals with the resummation in the 1 + 1 jet 
limit, of the distribution of the jet-broadening with respect to the photon axis, B^e, as 
measured in the current hemisphere of the Breit frame of deep inelastic scattering and 
defined as 

In the above equation n refers to the photon axis (conventionally taken to be the z 
axis) in the Breit frame and Tic is the current hemisphere. The above definition is 
valid provided one imposes a certain minimum energy cut-off Sq > i^um for reasons of 
infrared safety, with Sc the energy in the current hemisphere. The choice of Sn^ should 
be sensible (a not too small fraction of Q) to avoid development of further significant 
logarithms involving this quantity.^ 

In the 1 -|- 1 jet limit, the broadening is small and in the perturbative expansion of 
the distribution each power of can be multiphed by up to two powers of In B^e- This 
leads to a very poorly convergent series and necessitates a calculation which sums the 
dominant terms of the perturbative series at all orders — a resummation. 

As we have already mentioned, techniques for the resummation of event shapes in 
e~^e~ are well established and more recently there have been extensions to DIS jjlT 



to non-global variables |12| and to multi-jet configurations |15|,pi 

In general for a resummable variable V one can show an exponentiation of the large 
logarithms, which means that the suitably normalised cross section for the variable to 



^One could also replace the limit on £c with a limit on X^ieHc 1^*1 — ^^^^^ would reduce the 
sensitivity to non-perturbative effects associated with hadron masses [Q. A procedure allowing one 
to avoid the cut altogether would be to normalise to the photon virtuality Q. In contrast to the thrust 
case lll|, there would be no difference in the ensuing resummation. 
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be smaller than some value V can be written in the following form: 



a{V) 



(^ + Y1 '^""s j e^^i("^^)+^2(°=^)+- + D{V) , (1.2) 



where L = \nl/V. The distribution is obtained by differentiating this expression with 
respect to V. 

The leading (or double) logs are those contained in Lgi{asL) and the next-to-leading 
(or single) logs are those contained in (72 (og-^)- Further subleading sets of terms would be 
contained in functions ag~^5'„(asL). There is also a remainder function D{V) containing 
terms which go to zero for V 0. For problems with initial state hadrons (DIS, pp) 
the Cn and D, as well as g2, gs, ■ ■ ■ are generally x dependent and involve convolutions 
with parton distribution functions (which in eq. (|1.2| ) have not been explicitly shown). 

The summation of leading and NL logarithms is usually sufficient to control the 
normalisation of the distribution down to the region of the peak of the distribution, 
agL ~ 1. The contribution from further subleading terms, such as assy's («s-^); will be 
of order ctg in that region and so formally a 'small' correction. 

The variable considered here is however unusual in that if one follows the standard 
procedure and keeps just the leading and next-to-leading logs, calculated in section ||, 
one obtains an answer which diverges for agL of order 1. Such an occurrence is not 



unknown in problems involving initial state hadrons. In Ref. it was observed for 
the case of Drell-Yan lepton pair production near threshold, that a factorial divergence 
occurs (shown to be unrelated to the expected renormalon behaviour) in the coefficients 
of the resummed formula if one performs a naive resummation till NLL accuracy. 

Additionally, problems closely related to those discussed here have been observed 
before in certain approaches for calculating the transverse momentum distribution of 
a Drell-Yan pair |T8|, p!9|] , but to our knowledge the DIS broadening represents the first 
time these difficulties have come up in the context of an event shape — i.e. a variable 
with direct sensitivity to hadron emission. 

Essentially the divergence stems from the fact that there is a dependence of the 
observable on the transverse momentum recoil. In close analogy with the case of the 
Drell-Yan pt distribution, the net vector sum of the recoil can go to zero in two ways: 
either by a veto on emissions, or by vector sum of several emissions adding up to 
zero ||2^ . The exponentiated form for the answer is suitable for taking into account the 
first effect but not the second, and breaks down (with a divergence) when the 'easiest' 
(most likely) way of producing a low transverse momentum recoil is through the second 
mechanism. 

There are however important differences between the Drell-Yan case and B^e- The 
Drell-Yan transverse momentum is sensitive to emissions exclusively through their re- 
coil. The broadening however is sensitive to emissions in the target hemisphere only 
through recoil, but to emissions in the current hemisphere through both recoil and direct 
'measurement' of the emissions. This means that standard solutions for the Drell-Yan 
problem [|I^,|2D| are not so easily applied to B^e- 
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Accordingly in section |^ we develop a new technique for supplementing the answer 
so as to avoid the divergence. Essentially we redefine our accuracy criterion in terms 
of the relative impact of a given contribution on the final answer rather than in terms 
of a formal counting of logarithms. Technically this requires the expansion of certain 
integrands to be carried out about a point closer to their saddle-point than is needed in 
the standard approach. We show that after the application of this method remaining 
uncertainties are pure (non-divergent) subleading terms. 

For actual phenomenological applications, to be able to study the event shape dis- 
tribution over its full range, it is necessary to match the resummed calculation to 
fixed-order results. There exist well established techniques in e^e~, however for a vari- 
ety of technical reasons they cannot be directly applied to the DIS case. So in section § 
we examine the modifications that are necessary as well as elaborating on a matching 
scheme proposed in (here named M-matching) and introducing a new scheme which 
we call M2 matching. 

A final element in the prediction of event shape distributions is the non-perturbative 
correction. For most variables this is quite straightforward, being a simple shift of the 



distribution over most of its range |21,22|. However for the broadening, as is calculated 



in section ^, in close analogy with the e+e~ case I^S], the effect of the non-perturbative 
correction is also to squeeze the distribution. One interesting consequence of this is 
that the power correction to the mean broadening acquires an x-dependent component, 
which has been noticed experimentally by the ZEUS collaboration [[10|] . 

Given all these ingredients we are therefore able for the first time, in section ^, to 
show a comparison of a resummed, matched and power-corrected distribution to DIS 
event shape data . In that section we also show how the resummed results compare to 
the exact fixed order calculation and examine the effect of the standard and improved 
resummations. Forthcoming work will give a more detailed analysis of the data. 



for a range of observables including the broadening. 

In addition to the contents of the body of this paper described above, there are sev- 
eral appendices containing details of the working of the various sections. One appendix 
which we draw particular attention to is appendix ^ which deals with the evolution of 
parton distribution functions (PDF). Though PDF evolution is not the subject of this 
paper, it turns out that for the phenomenological implementation of our formulae there 
are certain advantages {e.g. freedom in the choice of the value of for the evolution) 
to using one's own PDF evolution rather than that embodied in standard PDF global 
fits such as those from CTEQ, GRV or MRST p^-p6|. Appendix |F| discusses these 



advantages in detail, presents the algorithms used, and also shows some discrepancies 
(though fortunately in contexts of limited phenomenological importance) that we have 
found in the evolutions embodied in the CTEQ5 and MRST99 distributions. 



2 Derivation 

Many aspects of the resummation of event shapes have become standard in past years. 
Accordingly in this derivation we will shall be quite concise, referring the reader to the 
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literature 



for a more detailed discussion of certain subtleties. 




We start off by writing the momenta of radiated partons (gluons and/or quark- 
antiquark pairs) in terms of Sudakov (light-cone) variables as (Fig. 1) 

ki = aiP + (3iP' + kti, aA = kl/Q\ (2.1) 

where P and P' are light-like vectors along the incoming parton and current directions, 
respectively, in the Breit frame of reference: 

P = xPp, P' = xPj, + q, 2{P ■ P') = -g2 = q2 ^3.2) 

where Pp is the incoming proton momentum and x = Q'^ /2{Pp-q) is the Bjorken variable. 
Thus in the Breit frame we can write P = ^Q{1, 0, 0, —1) and P' = ^Q{1, 0, 0, 1), taking 
the current direction as the 2;-axis. Particles in the current hemisphere Tic have (3i > ai 
while those in the proton remnant hemisphere Tip, have > 

In line with the procedure adopted in |Tl|, we write the following expression for 
the cross section (strictly speaking the contribution to F2 from an incoming quark of 
unit charge). It is given in terms of iV, the moment variable conjugate to Bjorken-x 
for the case of emission of m gluons in the remnant hemisphere and n in the current 
hemisphere, where all the gluons have kt <^ Q, 

(2.3) 



with the coupling defined in the Bremsstrahlung scheme |2^ and 



PR, = I dz, \^ QiQa, - ku), (2.4a) 
PC, = / dzi e(Q/5, - kt^) . (2.4b) 

1 Zi 
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Because of the coUinear divergence on the incoming (proton) leg, we have had to intro- 
duce a factorisation scale, Qq. This scale serves both for the parton distribution and as 
the lower limit on transverse momentum for emissions in TYr. The virtual corrections 
are given by expressions which are similar except for the absence of the factor zf in 

In the situation where there is an ordering > kt2 > ■ ■ ■ hi, the and Pj can be 
written as follows 

a. = ^^, Pj = ^^. (2.5) 

Zi . . . Zi Zi . . . Zj 

For differently ordered cases one should simply permute the indices appropriately. In the 
limit of soft emissions the Zi can all be approximated by 1. As is explained in to our 
accuracy it is actually possible to do this even when there are hard collinear emissions. 
So the factor Q{Qai — ku) can be replaced by Q{Q{1 — Zi) — ku) and analogously for 

0(g/5,-fcy). 



2.1 B 



zE 



Here we work out the contribution to F2 from events with a broadening smaller than 
BzE- We shall examine configurations consisting solely of soft and/or collinear emitted 
gluons. In this limit the difference between 2 Ylii^Hc 1^* I ^^"^ ^ small and given that 
the broadening is also small we can replace 2 ^ .^^^ \pi\ Q introducing an error of 
O {k^/Q"^) which is negligible. In order for the broadening to be smaller than Bze, one 
then obtains the following condition on the emitted momenta: 



\ He 



(2.6) 



The real emission part of the contribution to F2 (for an incoming quark of unit charge) 
is then given by 



'^n{Bze) = ^ CTmn© j QBzE — ^ hi 

-He 



(2.7) 



In order to carry out a resummation we need not only (Jmn in factorised form, as given 
above, but also the B function. This can be obtained with the aid of a couple of integral 
transforms: 



e 



du d^bd^pt 



n 



-ukti/Q 



n 



(2.8) 



where we have explicitly introduced the current quark transverse momentum pt 
- J2ienR,nc Carrying out the pt integration gives [Q 

du Q'^d% V 



e 



2mv 271 [u^ + h^Q'^fn 



n 



-vkti/Q 



n 



Ab-kt. 



2.9) 
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We then define the double transform ^7v('^, b) of our cross section: 



2TTiiy 2tt 



-^Ar(l^, b). 



(2.10) 



Writing the sums over m and n as exponentials one obtains the following all-orders 
expression for <;n{^, b): 



^Ar(z/, 6) = exp 




In 



It: 



1 _ 

^ Q 



dz^^^iz^'e''-'^-! 
1 — z 



1-- 



dz 



1 - z 



>.fcV.WQ_i] \q^{Ql). (2.11) 



We have simplified the phase-space restrictions in eq. (|2.4|) , making the approximation 
that zi . . . Zi ~ 1, as explained in above. The terms (—1) account for the virtual 
corrections We then write 



(!n{QI 



exp < - 



nkf 



2n 



dz^jt^ {z'' - 1) \ qM') , (2.12) 



where, as elsewhere so far, we keep only the piece relating to gluon emission from a 
quark. We have exploited the fact that whenever an integrand is finite for z —>■ 1 we 
are allowed to ignore the (sub leading) difference between 1 — kt/Q and 1 in the upper 
limit on z. It is then convenient to rearrange the iV-dependence so as to isolate the 
■Jqq,N anomalous dimension: 



dz- 



Z L 



1) 



..N 



1) 



Ab.kf 



lqq,N + 



dz 



z 

(2.13) 



Since the A^-dependence is now completely separated from the 1/(1 — z) soft divergence 
it is straightforward to replace 7qg,Ar with the full anomalous dimension matrix, 7^. 
Accordingly our complete answer in Mellin transform space is: 



?Ar(z/, 6) = Cq^nG 



(2.14) 



where Co,7v is a matrix of zeroth order coefficient functions in space (see |jTT|]), q is 
a vector of parton distributions and one has the following expressions for radiators TZc 
and T^r: 



7^R(&) 

7ec(z/, b) 



dkf Cpasikf) / 
dkl Cpa^ik'i) I ^ , ^-vktlQ 



Q , (l + 2;2^ 
dz- 



271 



{Mbk, 



t e 



1) 



dz 



1 - z 

1 + z^ 
1 - z 



(2.15a) 
(2.15b) 
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To NLL accuracy these integrals can be evaluated by making the following replacements 

EM 



{Uhh) e-^^^ - 1) - I fc, - ^7W^ I ' (2.16a) 
{Joibkt) -l)^Q{kt- . (2.16b) 



b 

Carrying out the integrations over z and z one obtains 



7^c(.,6).i?(^±V^^±^l, (2.17a) 
7^R(6) ^ i?, ( ^ ) , (2.17b) 



2 

where we have introduced u = ve^^ and h = be'^^, and 

Explicit expressions for R{u) and -R^(m), to leading and next-to- leading accuracy, are 
given in appendix ^ 

Our answer in u space is then given by an integral over b: 

= C„,„y MS^^^j^pJ^e { ' ) <'>,„(Q'). (2.19) 



To evaluate this integral, the procedure that has being adopted previously 0,[Tl|], has 
been to expand the functions R{u), eqs. ( p.l7[ ), as 



R |^^±V:^±W j = ^(2z>) + i?' In l±^y^ + O {R" In^ y) + O {R', In y) , 

(2.20a) 

^"(?) =^7(^) + ^;in| + C?(^^;in'2/)+(^K2lni/) , (2.20b) 

where we have introduced y = bQ/v. We have defined R! as been the pure single- 
logarithmic piece of vdyR{v)\ accordingly R!^ = R! . The part of vdyR{v) containing 
terms a" In'^"^ z/ is referred to as R'2 (and analogously for R^). In what follows imme- 
diately below it can be neglected because it leads only to NNLL corrections. The same 
criterion also allows us to throw away the terms containing R" . 
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So we can now write our expression for the v-space resummed cross section: 

?jv(^) = Co,N e-^-(^)-^(2'')A(2i?') qM') , (2-21) 

with 



-R' 



-R' 



We are not aware of a closed form for A. Its expansion for small R' is 

A(2i?') = 1 + l!L!zJ^l!^(2i?')' + O ( R'') . (2.23) 
96 V / 

The final step of the calculation is to perform the inverse Mellin transform with 
respect to u: 

^n{B)= I ^e^VM- (2.24) 
It can be shown in a variety of ways (see for example [^]) that this gives 

^NiB) ^ — ^ / , rc;^ (^) . (2.25) 

r(l - z/9^1n^A,(z/)|^=i/B) \BJ 

Accordingly after throwing away all NNLL (and yet higher order) terms we obtain 

Its fixed order expansion is given in terms of the coefficients Gnm of (as/2vr)" In'" 1/B 
in the exponent, which are listed in table |l|. Note that in contrast with the convention 



adopted in |Tl[] the change in scale of the parton distribution has not been written 
explicitly, but rather left implicit through the action of R^ on q^r. Accordingly the 
Gmm are not pure numbers but rather operators in flavour space. 

We note that the answer has been checked against an analytical first-order calcula- 
tion of the dominant terms at small B, given in appendix It has also been tested 
(strictly the form ( [4.2|) , which includes the O (as) constant term) against fixed order 
results from DISASTER-|--|- and we find good agreement for all terms that are 
intended to be under control, namely cxsL^ with < n < 2 and a^L" with 2 < n < 4. 

Finally, we observe that while we have chosen to resum logarithms of 1/-B, we could 
just as easily considered a resummation of say logarithms of 1/XB. With an appropriate 
change of the g2 function, given in appendix y one then obtains an equivalent answer to 
NLL order, but with different effective subleading dependence on B. Such a rescaling of 
the argument of the logarithm is therefore a way of testing the sensitivity of the answer 
to uncontrolled subleading effects. 
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Gl2 




Gil 


(6 -41n2)CF - 27^ 


G23 


3 '^F 


G22 


[{12 - 16 In 2)7r/3o -47^-^^ + 24 In^ 2) Cp) Cp - 4vr/3o7jv 



Table 1: Fixed order coefficients of the resummation for B^p., at order ag and 



2.2 Problems 

Though eq. (|2.26|) is correct to NLL order, it turns out that even in the region where 
the expansion is formally valid, i?' ~ 1, there are problems. For R' > 2 the integral 
( |2.22| ) diverges in the y —>■ region. One can examine in more detail what is happening 
by including the R" term of our expansion around bQ = u in eq. ( p.20b| ). One finds 
that for R' close to 2, the modified form of the integral (|2.22| ), gives 



i. e. the formally subleading R" term (as well as yet higher-order terms of the expansion) 
becomes enhanced and can no longer be neglected. 

The breakdown of our expansion around b = u/Q is associated with two physical 
facts: ffistly, half of the double logarithmic contribution comes solely from recoil of the 
current quark with respect to emissions in TYr; secondly the recoil transverse momentum 
of the current quark can be zero even if it has radiated gluons — it suffices that the 
vector sum of the emitted transverse momenta be zero. 

It has been known for some time |]2D| that there are two competing mechanisms 
for obtaining a small recoil transverse momentum. One can restrict the transverse mo- 
mentum of all emissions in TYr — locally the probability of getting a small transverse 
momentum from this mechanism scales as (there is another factor pf coming from 
a restriction on emissions in Tic, however this factor persists independently of any dis- 
cussion of recoil because it is also generated by directly observed gluons). Alternatively 
one can have a small recoil transverse momentum due to the cancellation of larger 
emitted transverse momenta, and the corresponding probability scales as p"^. 

In the region where R' < 2 the easiest (least suppressed) way of restricting the pt is to 
restrict emitted transverse momenta. Accordingly one sees a probability proportional 
to a Sudakov form factor. However for R! > 2 this Sudakov form factor associated 
with Tipj, emissions gets frozen (at its R' = 2 value) and the alternative mechanism of 
suppression takes over. 

The divergence that we see is associated with this transition. It arises because we are 
trying to use a formula, eq. ( p.21 ), with the same double logarithmic Sudakov structure 



above and below R' = 2, and a single logarithmic factor A{2R') intended to account 
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for the effects of multiple emission. However there is no way for a single logarithmic 
function to cancel the effect of a double-logarithmic Sudakov form factor and bring 
about the "freeze-out" discussed above — so at the point where this is supposed to 
happen, the simple approach of eqs. ( |2.20| - |2.2^ ) breaks down, giving a divergence. 



Mathematically, what this corresponds to is that for very small B the integral ( ^.191 ) 



is dominated by values of <C z/ (related to individual emitted transverse momenta 
in TYr being much larger than BQ), so that an expansion around bQ ^ i> is bound to 
fail. 

Related problems have been seen before, in certain approaches for the calculation of 
the transverse momentum distribution of a Drell-Yan pair [|18],|19|. They turn out to be 
a general feature of observables for which the contributions from different emissions can 
cancel. Other examples of such variables are the oblateness (the difference between the 
thrust major and minor) [^] and the difference between jet masses in e~^e~ . Strictly 



speaking even the thrust Tze, resummed in |]Tl|, suffers from this problem — however 
there, for actual values of jSo, the divergence turns out to be to the left of the Landau 
pole and so can be ignored. 



3 Beyond the divergence 

For most phenomenological purposes it turns out that the divergence at fZ' = 2 does 
not cause any practical problems. This is because it is considerably to the left of 
the maximum of the distribution {R' = 1/2), in a region where the distribution is 
strongly suppressed by the Sudakov form factor, and where there are in any case large 
uncontrolled non-perturbative corrections. 

However one can envisage cases where the divergence may cause problems (for ex- 
ample when using a non-perturbative shape function so as to extend the distribution 



down to zero B [^,^), especially at low Q values, where it is more pronounced. Fur- 
thermore it is in a region which should formally be under control. So we feel that it is 
worth dedicating some effort to improving the answer in this region. Additionally the 
techniques that we develop may be of use for other variables where similar divergences 
occur closer to the phenomenologically relevant region. 

Before entering into the details of the method it is perhaps worth commenting on 
criteria for including subleading terms. In the standard approach it is usual to keep the 
minimal set of terms — i.e. just the leading and next-to-leading logs, and throw away 
anything which contributes beyond this accuracy. This is analogous to the philosophy 
in a fixed-order calculation, where one keeps only the orders one knows and sets higher 
order terms to zero. Accordingly for example, whenever we have an expression involving 
OlR we keep only the dominant (single- log) terms R', eq. (|A.7|) , because other terms 
would lead to next-to-next-to-leading contributions. 

For a 'normal' observable, in the region a^L ~ 1, the inclusion of NLL terms is 
required to guarantee that the relative error on the answer is of order ag (associated with 
the NNL aggsidsL) term in the exponent). For the broadening a natural extension of 
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our accuracy criterion, is to keep not just a particular set of logarithms, but additionally 
all terms whose relative contribution in the region agL ~ 1 is larger than O (as), even if 
they are formally N'^LL, with n >2. In general we will try to take prescriptions which 
are as close to the original formulation as possible. We refer to these as 'minimal' 
prescriptions. 



3.1 Improved resummation 

In the standard approach the divergence arises only when taking the inverse Fourier 
and Mellin transforms of our answers. Accordingly up to eq. (p.l9|) the resummation 
procedure remains unmodified. It is in performing the inverse transforms that the 
method needs to be improved. Roughly speaking our approach will involve expanding 
TZnih) around a point 60 close to the saddle point of the integral ( |2.19| ), rather than 
around h = Q/v. It will also be necessary to keep a larger number of terms in the 
expansion. 

So we start off by finding the point feo of eq. ( p.l9|) around which to expand TZr. To 



within our accuracy, as is shown in appendix |D|, it is necessary for to be close to the 
saddle point of the integral, but it is allowed to differ from the true saddle point by a 
factor of order 1. This enables us to choose &o such that in the limit of small R' we 
have h^Q = u, as in the standard resummation. 

A convenient way of doing this is to define 60 = Vo^/Q through the saddle point, 
I/O, of the following integral 

^m = - [ — expi-iL + \ny)R,ia,iL + \ny))) , (3.1) 
IT J y cosh 2 my 

with L = Inu. This corresponds to finding the solution of the following equation, where 
we have defined £ = In 7/0 : 

2tanh2^ + i?'(as(L + ^)) = 0, (3.2) 

or equivalently 

2 tanh 2i+^ = . (3.3) 

This prescription for choosing 60 is minimal not just because for small R' it gives 
i>f)Q = ^1 but also because all unnecessary subleading terms have thrown away. The 
solution to (|3.3| ) has the following expansion 



£ = -CpasL - CfL{A'k(3qL - Cf)^^ + ^ [alL^) , (3.4) 
where we have defined as = as/27r. We also define 

L = L-i, R' = R'{asL), R" = R"{asL). (3.5) 
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We then expand 7^_r(6) around 60 = , 



Tl„(b) = P^{ve') + 7,fl' + ( In ^ - C ) (fl' + i?;^ + 7,fl") 



+ if,„!!«_,yfl" + if,„^_A'ff", (3,6) 
2 V 2z/ / 6 V 2z/ / ^ ' 



Comparing to the expansion ( |2.20b| ), the additional R!' term on the second hne is 
required in order to control the answer to within a factor of O (1), while the ^275 7e-R" 
and R'" give contributions of relative order ^/a^. This too is discussed in detail in 
appendix ^ For IZci^, h) one expands as before, 



nc{y,h) = R{2u) 



(3.7) 



since the series remains well-behaved in the limit 6^0. 
We then write our z/-space answer as 

where 

-R' 

^ 1 r 11^ / T _|_ , /I _|_ \ 



(3J 



dy (1 + VTTF^ 



1 - 7 1/(1 + y^f'^ 



4 



g-(ln|-^)ij'-i(lnf-£)2i?" ^ 



X 



1-ifln^ 
6 V 2 



R!" - (in I - (i?;^ + l.B!') 



. (3.9) 



In the expansion of the exponent, terms whose relative contribution is of order ^Joi^ 
have been kept only at first order, whereas other terms must be kept at all orders. 
The factor 1/(1 — vr^agCF/S) is included so as to ensure that A is free of any O {a^) 
contribution in analogy with the standard resummation. Writing its expansion as 



A = l + ^^A^„<i^" 



m=2 n=0 



the order a1 terms are given by 
3 



A. 



20 



10 

4 



Cf'a^TT^ + C(3)[67jv + (127e - ^)Cf - 87r/3o]C^ , 



A21 = -Cf (-2 TT^/^o + 9 Cf C(3) + 6 Cf X12) + 



+ C^(27jv + (47e - ^)Cf) ( -1 + ^ - 2 In^ 2 



A22 = 7TC^(77r2-36 In' 2 -6) 



(3.10) 



(3.11a) 

(3.11b) 
(3.11c) 
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with 



X,, = I , 1n f i±V^ ] y ^ 1.945031318. (3.12) 

J y (1 + 1/2)3/2 1^ 4 y 2 ^ ^ 

We note that eqs. ( ^.91 ) and ( ^.11| ) involve the anomalous dimension matrix 7^, through 
its presence in the term ^1^2- 

It is useful to study the behaviour of the various factors of ( |3.8| ) in the two regimes, 
R' < 2 and R' > 2. For small L, the series expansion of i, eq. ( p.4|) , shows that i 
remains small compared to L. Indeed R-y{i>e^) differs from R-^{v) by a term of order 
ctgL^ (and higher) which is exactly compensated by the O {ajL"^) difference between A 
and A, cf. eqs. ( |m^ ) and (g]^. 



If one increases L, then one finds that a transition takes place around R' ~ 2. 
Beyond this point i starts to vary much more rapidly, going roughly as 

i^-L+ +0{1) , {R'>2) (3.13) 
(Cf + 27r/5o)as 

Accordingly R^i^ue^) stops varying for R' > 2; on the other hand A starts varying 
rapidly, going as e"^^. 

In section we mentioned the presence of two competing mechanisms for the 
current quark to have a small transverse momentum. The transition in the behaviour of 
i that we have just discussed corresponds precisely to the transition from the mechanism 
of suppression of radiation (associated with a Sudakov form factor), to that of arranging 
for the vector sum of the emitted momenta to be small (the probability of which scales 
as ) . 

It is important to understand this transition, and in particular the L dependence 



of the various factors in (|3.8|) in order to perform the inverse Mellin transform of our 
result. As in section ^ we shall make use of eq. ( [^.251 ). To calculate the argument 
of the F-function we need to know which factors in (|3.8|) vary rapidly. Since we aim 
to control our answer to a relative accuracy of O (as) this means that (logs of) terms 
whose derivatives are of order a^L, or of order 1, must be kept, while terms for the 
which the derivative is of order can be neglected. 

So whereas in section |^ this meant that we could neglect A, now, since A can vary 
rapidly it needs to be taken into account in i/c^j, ln^Ar(i/). Nevertheless, in line with the 
approach of neglecting non-essential subleading terms, we shall not take the derivative of 
the whole function ln^Ar(i/), but only of those terms which are essential. One possibility 
would be to define 

KR = -^iLRi{a,L)~lnA) (3.14) 

This however presents some technical difficulties, because A involves the anomalous 
dimension matrix. These difficulties could perhaps be surmounted, however a simpler 
solution is to observe that while R'^2 contributes pieces of relative order ^/a^, they 
vary significantly only over a region of AL ~ l/^/a^ (see appendix |^). Accordingly in 
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Sl In A they can at most contribute an amount of order as and so can be neglected. So 
to avoid the comphcations associated with the anomalous dimension matrix, in ( p. 141) 
we could replace A with the following 'simpler' quantity: 



A„ = ] [ ( 1 + ^1+^^ j g-(ln|~^)i?'-i(ln|-^)2ff 



-R' 



l-i^a.CpJ y (1 + 2/2)3/2 4 



X 



X 



. (3.15) 



which differs from from A only through the replacement of B!^2 -^2 i'^ ^^e second line. 
Correspondingly its fixed order expansion differs from that of A by the absence of the 
7jv terms in ( p.lla| ) and (|3.11b| ). The subscript n indicates that this is a non- minimal 
choice for A. 

We can also make a more extreme choice, throwing away all terms which do not con- 
tribute significantly to the derivative. One possibility makes use of the same integrand 
as was used to determine the saddle-point of and gives 

where the subscript 'm' stands for 'minimal'. The terms of its fixed order expansion 
that will be needed are 

A...^-.%C., A„.^(^-4)CJ. ,3.17) 

Our final answer for the integrated broadening distribution will therefore be given by 

where B = e~^B and 

Ks,x = ^ {LRi{a,L) - In Ax) , (3.19) 

with X = s or m. Here, we have written the scale of the parton distribution explicitly, 
to emphasise that it is now BQ rather than BQ. 

For the purposes of matching to fixed order we need to know the expansion of ( p. 181 ) 
to 0{ag). One finds that it differs from that of ( p.26|) by the following additional 
(subleading) terms: 

[A20 - Ax2i7e + (A21 + (277V - SCf)Cf - (2Ax22 + 4C|)7e) L] . (3.20) 
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4 Matching 



In order to extend the range of validity of the predictions for event-shape distributions, 
various procedures have been developed in e~^e~ |1| for supplementing resummed distri- 
butions with the information from the first and second fixed order distributions. This 
is referred to as matching. It is not possible to merely carry over the matching schemes 
developed in e'^e'^ to DIS without addressing certain technical comphcations that arise 
in the DIS context, which shall become evident here. Additionally we propose new 
schemes that are suitable for the purposes of matching our resummations to fixed order 
estimates. The discussion that follows will be kept fairly general since we intend to use 
the various schemes and ideas introduced here not just for the jet broadening but for 



other DIS variables that have been resummed thus far [|rT],[T^. For this reason we refer 
to all distributions as being a function of a generic variable V rather than specifically 
of the broadening B. 

We begin by examining the form of our resummed cross section. This has the 
following structure in moment (A^) space 

^^r,7v(Q^ V) = [Co,iv + a,Ci,^] exp[Lg^{asL) + g^^asL)] qM^), (4.1) 

where the subscript denotes Mellin transformed quantities and we have introduced 
flavour vectors and matrices in 2n/ + 1 dimensions for the constants, the parton dis- 
tributions (see for details) and the functions gi and g2- The operator structure for 
Qi is trivial (diagonal matrices with the scalar function gi, computed as in appendix ^ 
for the broadening case, for the quark entries and zero for the gluon entry) and only 
needed for dimensional consistency. The operator structure of g2 is non-trivial due 
to the presence of the anomalous dimension operator. Note that this result is for the 
un-normalised cross section. Since a normalisation in A^ space does not correspond to 
a normalised x space result, we only normalise our result after translating to x space. 

The normalised resummed result in x space then reads as 



a(x V^O"^) 

ar{x, Q\ Y) = S(x, Q^ V) + a, exp[Lg, + g;]C, ® ^^^^^ ^ (4.2) 



with the form factor S defined by 

E(x, Q',V) = eMLgi + 9l f^2r^? ' ^^'^^ 
where we have introduced the singlet distribution 

q{x.Q^)= Yl e?fe(x) + g,(x)], (4.4) 

j=u,d,s,... 

and used it to normalise the x space result. Note that the piece of g2 corresponding 
to DGLAP evolution (that involving the anomalous dimension matrix) has been used 
to change the scale of the parton distributions to V^Q^ leaving behind the g2 function. 
All quantities in ([4.3|) are now scalars rather than operators since in writing the above 
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we have multiplied out the matrices involved. The result for Ci is available in appendix 
0. The index n can have the values (for variables like the jet mass, C parameter or 
thrust with respect to thrust axis), 1 (for the thrust with respect to the photon axis) 
or 2 (for the broadening). 

A point that we wish to draw attention to is that although the form factor contains 
a parton distribution evaluated at scale V^Q^, we can ignore this change of scale in the 
second {O (a^)) term of eq. ( ^4.2| ) (and use for the scale of the parton distribution) 
since it leads to subleading terms starting at the a^lnl/V^ level. However such a term 
will of course be relevant in the matching to NLO. Hence if one chooses not to keep 
the scale as V^Q^ in the convolution involving Ci we will have to modify the matching 
piece accordingly. 

Now we are ready to match the resummed result to the fixed order result returned 
by the NLO DIS Monte Carlo programs. Let us denote these exact results by cxe^'* and 
o"P^ for the and a1 Monte Carlo estimates. These are the result of a convolution 
with a specified structure function and so are returned in x space, and taken normalised 
to g(x, Q^). In order to perform the matching we essentially have to add the Monte 
Carlo and resummed results and remove the pieces which would be double counted. 
These pieces would be the O (ag) and O {a1) terms of the resummed result cxr^^ and 
CTr^'*, which can be obtained by expanding eq. ([4.2|). 

However note that the matched resummed cross section has to satisfy certain re- 
quirements. The most important property is that in the V — limit the cross section 
must vanish on physical grounds Accordingly the following matching formula is 
invalid 

aiV) =ar + as (^^^ - o^^^) + al (af - af^) , (4.5) 

('2') ('2) 

because for \^ — 0, the factor (de — (7r ) does not vanish, but rather grows as In 

In e^e~ the two main matching procedures are R and hiR matching [R in the 
original papers is the equivalent of a here). In R matching one determines (from the 
fixed-order distribution) the G21 and C2 coefficients for the distribution, and defines an 
improved resummation formula 

aR = {l + Cias + C2as)e^^i("-^^)+^^("=^)+^2^"'^ . (4.6) 

Then the equivalent of eq. ( |4.5|) with replaced by an does indeed vanish in the limit 

V 0. This procedure is feasible in e~^e~ because C2 and G21 are simple constants, 

(2) (2) 

and they can be evaluated by subtracting cxr from de in the very small V limit. 
However in DIS C2 and G21 both have x-dependence, and it is simply not feasible 
to extract numerically them with their x-dependence. One might think of extracting 
them (individually for each x, point) once the convolution has been done with the 
structure functions. However experience in e~^e~ shows that one needs to go to very low 
values of V, with vast statistics, in order to reliably extract these quantities. In DIS 
with DISASTER-I--I-, low values of V are often not accessible because of cutoff effects. 
Furthermore the Monte Carlo statistical errors tend to be an order of magnitude larger 



16 



than for a similar number of e^e events, and a similar number of events in DIS with 
DISASTER++ takes an order of magnitude more time to generate. 

The In R matching approach is more easily extended to DIS. The philosophy of In R 
matching is to carry out the matching in the logarithm of the cross section rather than 
in the cross section itself: 

av = exp [Ina, + a, ((Inae)^^) - (Ina,)^^)) + al ((Inae)^') - (Ina,)^'))] , (4.7) 

where (In cTe)^"-' is the O (a") part of IncXe. 

Strictly, taking the logarithm of the cross section is a delicate operation because of 
the operator structure, which enters at NLL level, in particular in the coefficients G„n, 
and one cannot for example use Ini? matching exactly in the form prescribed in . 
However using 

In ae,. = ln(l + a.a'^} + alaf)) (4.8) 
and expanding the logarithm in each case in powers of one can alternatively write 

which still retains the correct O (a^) expansion as well as the LL and NLL terms. Hence 
the above form is the one that should be used for Ini? type matching in DIS. Note that 
further subleading contributions are of course not exponentiated as operators. But since 
they are beyond our accuracy, we are entitled to mistreat them, as long they do not 
lead to some particularly pathological behaviour. 



Another matching scheme we considered in |Tl| was the following which we shall 
now call multiplicative, or M-matching, 

aiV) =ar+ [a, (a« - a^^) + al (af - af^) - 

0^1 (^?^ - <^^^) (L'Gu + LGn)] S(x, Q^ V) , (4.10) 

where the presence of the form factor, S, ensures that the whole cross section does go 
to zero for K ^ and we have used Qnm for the x space versions of the resummation 
coefficients listed in Table 1. Note that the Qu in the above result involves matrix 
products and convolutions in x space. For example for the jet broadening one gets from 
Table 1, by inspection 

Gii{x) = (6 - 41n2) - 2Co ® P ® -, (4.11) 

q 

with P being a matrix of leading order splitting functions. The x space versions of the 
Gn,n+i coefficients are the same as the N space numbers mentioned in Table 1 and we 
do not distinguish them notationally. 

We can also define M2 matching: 

aiV) = ar + as (a« - a^'^) + aj (af - a^^)) S(x, Q^ V) , (4.12) 
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which exploits the fact that the O (as) term does vanish in the V ^ hmit. (Assuming 
of course that one has the correct Ci). This has some similarity to R matching in that 
the O (ctg) piece of the remainder is not suppressed by a form factor. 

If there is also resummation off a gluon^ then we need to modify the matching 
somewhat. Strictly one would want some way in the fixed-order contribution of sepa- 
rating out contributions associated with the presence of just two gluons in the current 
hemisphere (or one gluon plus virtual corrections). However with the existing tools, 
DISENT and DISASTER++, this is not possible. Accordingly we arbitrarily choose to 
attribute the entire difference between exact fixed order and the expanded resummation 
to the part of the resummation that is off a quark leg. The formulae for M and M2 
matching remained unchanged, while that for InR matching becomes (where arq, (Jrg 
refer to the resummations off the current quark and gluon respectively) 

cry = arg + exp [Incr^g + as ((In (Jeg)*^^^ - {Inarqf^^) + al ((Inaeg)^^^ - (Ina,.,)^^^)] 

(4.13) 

where we have defined a^q = a^. — a^g. After explicitly taking the logarithms we obtain 

fi („W . ^2('^(2) m 1/(1) / (1) , ^(1) 

a{V) =arg + (Trqe " J+°4"^ -n'-^ j^'^ 

There are some other important requirements of the final matched result. One is 
that at the upper limit of the distribution V^ax the integrated cross section must go 
to its exact upper limit without any additional leftover tems O (ag). Where the fixed- 
order differential distribution goes to zero at the upper limit one must ensure that the 
matched-resummed one does the same. In order to obtain these properties, one should 
use modified matching formulae. The details can be found in appendix ^ 



5 Non-perturbative effects 

For most event-shape observables, in the Born limit the principal consequence of non- 
perturbative (NP) corrections is a uniform shift of the distribution by an amount of 
order Aqcd/<5 121])@-Q This is because the effect of low- momentum radiation on the 
observable is independent of the configuration of the hard momenta in the event. 

For the broadening the situation is more complicated because the effect of low- 
momentum emissions depends critically on the configuration of the hard momenta. 
There are actually several effects. One is that it is possible to neglect the recoil from the 
non-perturbative emission because after azimuthal averaging it is washed out relative 
to the recoil from perturbative radiation. Accordingly the only NP correction to the 

^We encounter this situation in DIS observables such as current jet-mass, C parameter and the 
thrust defined with respect to the actual thrust axis as well as for light jet masses and narrow jet 



broadenings in e+e |p3|,pj[ 

•^More precisely it is to convolute the perturbative distribution with a non-perturbative shape 



function of width Aqcd/Q |21|, however as long as this width is much smaller than the width of the 



PT distribution, this effectively reduces to a shift. 
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broadening comes directly from the transverse momentum of low-momentum particles 
emitted into the current hemisphere. (This is similar to the effect which reduces the 
naively calculated power correction to the heavy-jet mass by a factor of two pTj^). 

But there is a second very important effect: soft gluons are emitted uniformly in ra- 
pidity with respect to the quark axis. However transverse momentum is measured with 
respect to the z axis. When gluons are emitted at very small angles to the quark axis 
the contribution from their transverse momenta is entirely cancelled by a longitudinal 
recoil of the quark. So only when gluons are emitted with an angle larger than Q^qi the 
angle between the quark and 2;-axes, do they contribute to the NP correction. 

This has been discussed in detail for the e+e~ broadening in [|^, and the techniques 
developed there can be carried through almost in their entirety. Accordingly we shall 
only outline the steps, following very closely the working of section 3.2 of |2^, and 
illustrating the relatively small differences. 



5.1 Power correction to the distribution 

After integrating over the rapidity of NP emissions, the NP contribution to the broad- 
ening can be written as 

p(^ln^ + r/o^ , r/o ~ -0.6137056 , (5.1) 

where is the transverse momentum (with respect to the z-axis) of the current quark, 
and V governs the overall magnitude of the power correction [Q: 

7^ ^ ^A^^|ao(/i/)-«s-/?o^(ln£ + ^ + l)| , as = «Ms(Q)> (5-2) 

with Q;o(/i/) a non-perturbative parameter (corresponding to the first moment of ctg up 
to some infrared scale [ij] which is postulated to be observable and process independent. 
The 'Milan' factor M. accounts for the non-inclusiveness of the observable fB^,B4|,[ 



The effect of the non-perturbative contribution on the distribution can be deter- 
mined by replacing 

in eq. (|2.8|). Since P is a small quantity, we are allowed to expand the exponential and 
we write 

^N{y) = + yVfN{y) + O {V') , (5.4) 

where the non-perturbative information function /Ar(i^) is given by 

Jo (i + r)"*/^ 



19 



We make use of the results (|2.17] ), expand in powers of i?', and introduce the differential 
representation of p3i, to obtain 



-a/2 



y^-a/2 
2, 



2 - 7e - ^70 + In 



+ 

1 + v^r+F 



2\ ,,-a 



2(1 + 

Evaluating the integral gives 

-_R(eTEe-9a)-i?(2eTEe-Si) 



(5.6) 



a=0 



/jv(iy) = Co,Are" 



where 

T(a) 
and 



^A(a) 



2 - 7e - + 5a + X(a) + J(a) - ^— ) ^^(g^) 



(5.7) 



a=0 



1 + Vi + y' 



-a/2 



A(a)yo (l + y^)3/^ 



-a/2 1+ y/l+y^ 

2(1 + 2/2) 



(5.8) 



l + ^TT^ 



-a/2 



A(a)(l + ?/2)3/2 



-a/2 



3a/2 



dy. (5.9) 



(l + ^)l+a_ 

Note the presence of the factor 8"/^ in the subtraction term, required for a proper 
regularisation of the integral in the a — > limit. 

After carrying the v integration and accounting for the effect of the da derivative, 
we obtain 



e7Ee-Ca)_/j(2eTEe-'^<') 



A(a) 



r(l + a) 



1 + -55(5, a) 



with 



55(5, a) = P ( In 4 - J(a) - J(a) - 2 + r/o + ^(1 + a) - ^(1) + ^^^^^ ^ 

£5 a, 



, (5.10) 

a=0 



. (5.11) 



Finally, having evaluated the derivatives with respect to a and 'undoing the expansion' 
with respect to powers of 1/Q we obtain 



(5.12) 



In order to better understand our answer we shall consider two important limits, 
i?' ^ and K 2. 
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Limit of R' — > Q. A useful cross check of the answer is to compare it with one's expec- 
tations for R' = 0. In this limit the perturbative broadening is determined entirely by a 
single emission. Half the time the emission will have been in the remnant hemisphere, 
implying B = Pt/Q and InQ/pt = In 1/5. The other half of the time the emission will 
have been in the current hemisphere and B = 2pt/Q and so In Q/pt = In 2/5. Taking 
the average one obtains 

6B{B, 0) = P (^In i + ^ In 2 + ?7o) • (5.13) 

Noting that /(O) = ln2 — 2 and J(0) = 0, we see that this agrees with the full result, 
eq. (1^. 



Limit oi R' — > 2. A second limit of interest is the point where A diverges, R' = 2. 
The techniques used here for calculating the power correction are analogous to those 
of section ^ for determining A. We know that in the case of A they break down 
pathologically for i?' ~ 2. 

In the case of 6B the method also breaks down, but in a less pathological fashion. 
The reason is that 6B involves ratios of divergent integrals, and so stays finite. In 
particular 

J(4) = 0, J^(4) = 127, (5.14) 

which leads us to the result that 

55(5,4) =p(^ln-i + ^ + r/o^. (5.15) 

In the case of A we devoted considerable effort to obtaining a correct answer in the region 
R' = 2, even though this is somewhat to the left of the peak. Our motivation for doing 
this was two fold. With more sophisticated models {e.g. involving shape functions) for 
non-perturbative effects, a perturbative understanding of that region may still be of 
interest. Furthermore the method may be generalisable to other observables for which 
the breakdown occurs much closer to the peak {e.g. for the difference between jet masses 
in e^e~). 

For the power correction however it is not clear that such an effort is warranted: (a) 
the techniques that are required are probably more complex than for A (even without 
a specific treatment of the i?' ~ 2 region, the working for the power correction is 
somewhat more complex than for the PT distribution); (b) in the region R! = 2 the 
simple approximation of a non-perturbative shift to the distribution is in any case 
thought to be a poor approximation; (c) any techniques developed would probably be 
useful only for the broadening. 

Accordingly for phenomenology we advocate the use of eq. ( |5.11|) . If one wishes to 



venture into the region around R! = 2 while bearing in mind that this is almost certainly 
not a safe endeavour, one can ensure that the distribution remains well-behaved by using 



21 



the following extrapolation for 5B beyond R' = 2: 



SB{B,a) =V (^n^ + ^ + T]o + Om95537{a-4:)^ , a>4, (5.16) 

where the first derivative around a = 4 has been determined numerically. 

Another region which deserves some discussion is that of large B. Normally in e~^e~ 
the expression analogous to ( |5.11| ) is used over the whole range of B. Of course for large 



B the expression is not valid, because it does not take into account non-perturbative 
effects from a base configuration with 3 or more hard partons. Nevertheless, phe- 
nomenologically this approach works rather well as long as one does not go beyond the 
3-jet region, where the distribution is in any case suppressed. In particular it is not 
unreasonable to expect non-perturbative effects to shift the distribution to the right 
even around the upper limit of the 3-jet region, though perhaps not by exactly the 
same amount. This is because the 3-jet upper limit is not the kinematical upper limit, 
and extra soft gluon radiation is free to increase the value of the event shape. (In the 
only case for which a calculation exists, the C-parameter at the 3-jet limit, C = 3/4, 
the power correction is found to be about half that in the 2-jet region [p^] ) 

In the case of B^e in DIS (as well as many other variables in DIS), the situation 
is different — the 2 + 1-jet upper limit {B = 1/2) is also the kinematical limit for 
any number of particles and extra soft radiation cannot increase the value of the event 
shape (though it can reduce it). Accordingly it makes no sense to shift the distribution 
by ( pID around B = 1/2. 

Of course we do not know what the right answer is. However one solution, which 
does at least preserve the property that the distribution should not extend beyond 
-Bmax = 1/2, is to replace 



B 



Pnp 



SB^5B=(l-[-—\ \6B, (5.17) 

where pnp is an arbitrary positive power (which we would expect to take of order 1). 
We refer to this procedure as a modified power correction, in analogy with the modified 
matching of appendix |^, though we note that the value of pnp used for the power 
correction does not have to be the same as p used in the matching. As for the p used in 
matching it should be varied so as to gauge the systematic errors associated with the 
arbitrariness of the procedure. 



5.2 Power correction to the mean 

We can use the above results to extract the power correction to (Bze)- It can be written 

as 

{6BzE)ix,Q') = J dB^^!-^^j^^6BiB,2R'iB)). (5.18) 
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As can easily be seen, the integral converges for \nl/B ~ ag . Since our aim is 
only to control pieces down to a relative order of ^/a^, it therefore is possible to make 
considerable simplifications to both T,{x,Q'^,B) and to 6B. We use 

E(PT)(a;, B)^{1 + aM^ix, Q^) + alL^G^z) e«-^^^' , L = \ri^, (5.19) 

ID 

where Qu is as defined in ( |4.11| ), and 

SB{B,2R'{B)) ^ 6B{B,0). (5.20) 

We then exploit the following relations (dropping the explicit x, labels, for compact- 
ness) 

dB^-— = l, (5.21) 



dB 

t/S^PT) 1 1 I ^ 1 g^^ 1 



— ln^ = ^J— f— -^^ + ^^ + 0(v^) , (5.22 



dB B 2V -Giaas 2 da 2 G' 



2 

12 



to obtain 



{5B,e) = -;\ —TT + 7 77^ ^ + ^o + P, 5.23 

where is to be evaluated at scale (or /i^). We recall that q is the vector of quark 
and gluon distributions, with q defined in ( [4. 41 ), and that P is the matrix of leading 
order splitting functions. An alternative form for ( p. 231) , which may be more practical 
to evaluate (but which introduces subleading corrections at O (og)), is the following 



= (ite ^ i - siTvrfW -w; + 'i« + o(vs:))v. (5,24) 

If one determines the derivative of the quark distributions numerically, then one should 
ensure that they are reasonably smooth in Q^, which as discussed in appendix is not 
always the case. 

We emphasise that the power correction acquires explicit x dependence, through 
the dependence on the scaling violations of the quark distributions. It is the first time 
that such a phenomenon is seen for the mean value of DIS event shape. It would be of 
interest to make a comparison with the results of the ZEUS collaboration, whose data 
seem to require non-negligible x-dependence [0 in the power correction. 

Finally we point out that the modification of the power correction described at 
the end of the previous subsection, ( p. 17] ), only affects the answer for the mean power 
correction, ( |5.23| ), at the level of terms O (y^), which are beyond our accuracy. 



6 Analysis of results 

Here we present numerical results based on the the calculations of the previous sections. 
All figures have been generated with as(M^) = 0.118 and where relevant ao{fii = 
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Figure 2: A comparison of fixed order, resummed and matched- resummed distributions, 
at two X, Q values for which measurements exist at HERA ||^. The standard resumma- 
tion has been used, and the matching is modified M-matching. No non-perturbative 
corrections have been included. 



2 GeV) = 0.5. Leading order parton evolution has been used, consistent with the 
philosophy of keeping only leading and next-to-leading logs in the resummation formula. 
The renormalisation and factorisation scales have been kept equal to Q. 



6.1 Resummation versus fixed order 

Figure ^ compares the resummed results (with and without matching) to the fixed order 
(ttg) prediction as calculated with DISASTER++ |2^. The two plots correspond to 



different Q values. One sees that in the small- -B region the resummation has a dramatic 
effect and that this effect is larger still at lower Q values. One also sees that the matched 
curves are essentially identical to the fixed-order results at large B, while at small B 
the effect of matching amounts to a small modification of the pure resummed results. 

Close to the upper limit of the distribution one sees a small secondary peak, most 
prominent at lower x and Q values. It seems that this structure may actually be an 
artifact of the interpolation of the fixed-order distribution, since there are arguments 
that suggest that close to the maximum the distribution has an integrable, — 5)^^/^, 
divergence: when there is a single particle in the current hemisphere, at an angle 6 
with respect to the photon axis, the broadening is B = ^ sinO. For 6 close to 7r/2 there 
is roughly a uniform distribution of 9 values, then this translates into a — S)"^/^ 
behaviour for the distribution of -B. In practice this divergence will almost certainly be 
smoothed out by soft gluon radiation, but that is beyond the scope of this paper. 

If we calculate the matched curves with DISENT (the program is much faster, 
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Figure 3: A comparison of the standard, 'improved' and 'minimal improved' kinds of 
resummations. The right hand plot shows the same curves as the left-hand plot but 
with a magnified abscissa. All curves use modified M-matching and are shown without 
non-perturbative corrections. 



but is known to give wrong subleading logarithms of B ||rT|) we find that the results 
are modified by an amount of the order of a percent (for the lower Q value). This 
small difference is a consequence of the matching that we have used: in the small-i? 
region the matching terms are multiplied by a Sudakov form factor, and therefore so 
are the discrepancies in DISENT. If the matching term were not multiplied by a form 
factor (as for example would be the case in R matching) then the discrepancy would be 
considerably larger, of the order of 6% in the peak region — however given the difficulty 
of carrying out R matching in DIS, this is unlikely to pose a problem. 

Finally we note that at the higher Q value shown there could be some non-negligible 
contribution from Z-boson (rather than photon) exchange. From the point of view of 
the resummation this makes no difference, but for the non-logarithmically enhanced 
parts of the fixed order calculation there could be an effect. Unfortunately of the fixed- 
order programs that are able to calculate the broadening distribution with reasonable 
accuracy, DISENT and DISASTER-I--I-, neither implements Z-boson exchange. 



6.2 Resummation variants 

Figure ^ shows the three different varieties of resummation that have been developed in 
this paper. For the standard resummation one sees some 'noise' in the left-hand plot. 
In the right-hand plot, which simply has a higher resolution on the B axis, one sees 
clearly the divergent structure associated with the derivative of the pole in A (which 
has been analytically continued to give it meaning beyond R! = 2). At lower Q values 
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Figure 4: Comparisons of the resummed results with HI data 0. Experimental statis- 
tical and systematic, added in quadrature (and averaged if asymmetric). Note that in 
the right-hand plot the x values are correlated with the Q values, going from x = 0.014 
at Q = 15 GeV, to a; = 0.18 at Q = 81.3 GeV, and that the distributions have been 
rescaled by a factor 10* with i running from to 5 respectively. 



the problem is more prominent, but remains confined to a region which though formally 
within perturbative reach, in practice is beyond the applicability of the formulas. 

We also show the two variants of the improved resummation. The minimal improved 
curve is actually very close to the standard resummation, while the non-minimal im- 
proved curve is somewhat different. This difference is indicative of the size of uncon- 
trolled subleading effects. If one wishes to use an improved resummation we recommend 
the 'minimal' variant, mainly because of its similarity to the pure LL plus NLL resum- 
mation. We note also that the non-minimal improved resummation can show some 
small instabilities (not visible here) for R' ^ oo, which arise in the derivative of the 
((1 + 1^1 + ?/2)/4)^^' factor in the integrand for A„. 

6.3 Comparison with data 

Figure ^ shows the comparison of our results to some data from HI [^]. The left 
hand plot shows a single x, Q value, illustrating the pure perturbative result and the 
prediction including the power correction (modified with pnp = 2). Because of the 
In 1/i? dependence, the effect of the power correction is not just to shift the prediction 
but also to squeeze it, in a manner similar to what is seen for the e~^e~ broadenings. 
Overall the description of the x = 0.11 data is reasonably good (we recall that we 
are using standard values of as(M^) = 0.118 and ao = 0.5). The apparently poor 
description of the two lowest B data points disappears in large part after integration 
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over the bin widths. 

The right-hand plot of figure § shows a comparison between the data and the power- 
corrected resummations (now integrated over bins) for a range of x, Q values. There is 
a worsening of the description at lower values of Q, but this is to some extent expected, 
due to the increased relevance of subleading corrections (both ETJ and higher 

orders in ag). 



Elsewhere we shall present a more detailed analysis together with that for a 
range of other variables, including a study of renormalisation and factorisation scale 
dependence, of matching scheme and rescaling dependence (c/. appendix 0), and of 
various other subleading effects. 



7 Conclusions 

For an accurate description of event shapes over the whole of the available phase space 
it is necessary to carry out a resummation of logarithmically enhanced terms and to 
supplement it with a non-perturbative correction. This paper deals with the broadening 
BzE and is one of a series addressing the resummation of a range of event shapes in the 
Breit-frame current hemisphere in DIS [pi]-[IB|]. 

In principle the broadening calculations involve a straightforward extension of pre- 
existing resummation mill and power correction techniques. In practice several 



subtleties arise. The standard prescription of keeping just leading and next-to-leading 
logarithms actually leads to a divergent answer. This divergence is associated with a 
change in regime: at a certain point as one goes to very small -B, the coefficient of the 
double logarithm is halved, because pure Sudakov suppression stops being the favoured 
mechanism for producing the small B values and is partially replaced by a cancellation 
of recoil between multiple emissions. 

In principle the divergence occurs inside the region which is supposed to be under 
control. Accordingly we have developed techniques for extending the resummation into 
the region with the divergence. It is then necessary to modify our 'accuracy criterion' 
— rather than being the correct treatment of LL and NLL terms, it becomes that for 
as In 1/5 ~ 1 the distribution should be under control up to (but not including) cor- 
rections of relative order ctg (for normal variables these two conditions are equivalent). 

In practice the divergence lies in a region where relative corrections of O (as) are 
multiplied by such large (purely numerical) coefficients that the distribution is in any 
case not well constrained. Accordingly phenomenology is restricted to a region of larger 
5, where it turns out to be sufficient to use the 'standard' resummation approach. We 
note however that the techniques developed here may be applicable to other variables 
such as the difference between jet masses in e^e~ where a similar divergence is likely 
to be present inside the phenomenologically relevant region. 

We have also made other developments which are necessary for practical phe- 
nomenology. In e^e~ there exist two standard techniques for matching the fixed order 
and resummed calculations. In DIS one of them (i?-matching) is awkward to apply be- 
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cause of difficulties in extracting the required fixed-order information from Monte Carlo 
programs such as DISENT and DISASTER-I--I-. Another procedure, logi? matching is 
applicable, but needs to be modified to address some subtleties that arise in DIS but 



not in e~^e . Furthermore we give details on the matching procedure proposed in |TT 
naming it M-matching, and also propose a new procedure, M2 matching. 

For the numerical implementation of our resummed formulae it was useful to develop 
a software tool for the evolution of parton distributions. A by product of this work was 
the discovery of bugs in the NLO evolution embodied in the MRST (at small x) and 
CTEQ (in the DIS scheme) parton distributions. In due course a stand-alone version 
of the evolution program will be made public. 



Another software development, which will be discussed in detail in |TH], is a tool 



which exploits factorisation to allow fixed-order Monte Carlo events at a single value of 
X, Q to be reused for other x, Q values. For DISASTER-I--I- in particular this allows a 
gain of an order of magnitude in speed and was essential for the generation of the fixed 
order distributions used throughout this paper. 

With the availability of these tools we have therefore for the ffist time been able 
to evaluate the resummed distribution of a DIS event shape including NLO matching 
and compare the results to data. We find reasonable agreement especially at higher x 
and Q values. More detailed analysis, including the study of other variables, will be 
presented in forthcoming work ||13|| . 

Finally we note that the programs to calculate the resummed broadening distribu- 
tion, including the matching and power correction, are available from the following web 
page: http : // cern . ch/ gsalam/disresum/. 
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A Formulae for the radiators 



We write the radiators eq. ( |!2.18D to NLL accuracy as R{l/B) = LRi{asL) + R2{asL) 
with 

R,ia,L) = (-2A - ln(l - 2A)) , (A.l) 
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and 



W)^^ln(l-2A)' ^.^2A + (l-2A)ln(l-2A) 



47r2/32 



1 - 2A 



CfI3i [ 2A + ln(l -2A) 1 , 2 



271(31 



1-2A 



ln^(l-2A) , (A.2) 



where A = ag^oL and L = In 1/5. We have defined the coefficients of the /3-function 
to be 



/3o 



IICa - luf 



17C\ - bCAUf - SCpUf 



and the constant relating the gluon Bremsstrahlung scheme to the ms to be 



(A.3) 



18 6 J ' 9^^' 
For the radiator including the anomalous dimension, Ry we have 

-R7I = Rl 5 



R-y2 — R2 ~ 



In 



ln(l - 2A) . 



We also define 



and 



d 



„w «s -3C^ as2K/?oC^^A-47r/5iC^Aln(l-2A) 
/?2(as, L) ^ j^ma.L)) = ^ + ^ .pl[l-2Xf 

The analogous results for the case with the anomalous dimensions are 



= -^2 + 



27r 1 - 2A 



Finally we shall need the leading parts of the second and third derivatives: 

2Cf as 



i?; = R" 



R'!; = R" 



TT (1-2A)2' 
TT (1-2A)3' 



(A.4) 



(A.5) 
(A.6) 



(A.7) 

(A.8) 

(A.9) 
(A.IO) 

(A.ll) 
(A.12) 
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B Fixed order result 



Here we calculate the first order coefficient function Ci, which appears for example in 
( |4.2| ). It is needed in a variety of contexts — for example in many (but not all) of the 
approaches to matching with fixed-order calculations, and also for carrying out com- 
parisons down to OfgL^ accuracy with fixed-order calculations. Calculating Ci proves 
more cumbersome in DIS than in e^e~ since instead of a pure number one obtains a 
function of the variable ^ = Q^/(2p.g), where p is the four-momentum of the incoming 
(as opposed to struck) parton. Additionally there are various contributions to be com- 
puted, i.e. transverse and longitudinal parts of graphs with an incoming quark as well 
as boson gluon fusion. It proves convenient to first compute the 0{as) result for events 
with broadening B^e > B because at leading order this quantity, doesn't get any 
virtual correction. Note that this is complementary to the final quantity we want, a'^\ 
which requires the selection of events with Bze < B. The following relation is therefore 
required (which follows from unitarity) 



(B.l) 



where a[ll is the total 0{as) cross-section for all events, except those with an empty 
current hemisphere which are excluded throughout, since it does not make much sense 
to define current jet observables in such a situation. 

The computation of (Tc^^ is relatively straightforward up to terms of order B, which 



we do not require, and after applying eq. ([B.l[ ) (at the level of the corresponding 
coefficient functions ) we obtain the following pieces relevant to the computation of ai^^ 



F2 quark contribution 



^,(s,0 = -5(1-0 



In- 



41n25-31n — - 41n21nS + 21n22 + 3-31n2 



1+e 



(1-0+ B i-e 

F2 gluon contribution 

j^^(B,o = -[e + ii-e)] 



ln(l-0\ , 6^2-2^-1 



1-e 



1 1-f 
4 In -- - 2 + 2 In —A 
B 4 



+ 



2(1-0- 



7r 



(B.2) 



^01-0- (B.3) 



Note that the longitudinal contributions are absent since they cancel in the two terms 
of eq. (BTT). Then we have aj^^ is given by 

2x 1 Cis 



ai'\x,B,Q') = ^^ I ^[CEJ',{B,0q{x/0+TfJ',{B,09{x/0], (B.4) 
q{x) 27r 4 



with 



[B.5) 
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In the above qj{x) and g{x) are quark and gluon distributions (for a quark with flavour 
index j and corresponding charge electric charge Cj). The colour factors are as usual 
Cp = 4/3, Tr = 1/2 and nj denotes the number of active flavours. From the above 
leading order result the constant (Ci) and logarithmic (Gn and G12) pieces can be 
easily read-off. In particular the logarithms are in agreement with those obtained at 
this order from the resummation method. The transpose of the matrix Ci in section 



IS 



Cf(x) 



[X) 

ix) 



{X 



with 



\ 



(B.6) 



(B.7a) 
(B.7b) 



Lastly we add that the above results are valid for the DIS factorisation scheme. To go 
the MS scheme one should add the standard ms scheme F2 coefficient functions to the 
above results exactly as in eqs. (A. 9) and (A. 10) in appendix A of Ref. |TI[ . 



C Rescaling the variable 

Usually the resummation of a variable V is defined in terms of Inl/V^; however we 
can just as well define it in terms of Inl/XV where X is some arbitrary number. If 
we do this then in the resummed formula we need to make the following replacements 
(note that the meaning of overlined symbols, such as L, is not in any way related to the 
meaning of barred symbols, L, introduced when discussing the improved resummation): 

L->L = ln^, (C.l) 

gi{asL) ^ gi{asL) , (C.2) 
g2{asL) g2{asL) + {gi{a^L) + a^Lg[{asL)) InX , (C.3) 

where bold symbols for those functions that are operators in x and flavour space. 
We also need to modify the constant term: 

Ci ^ Ci = Ci + 6-12 In^ X + Gn In X (C.4) 

In the fixed order expansion there are corresponding modifications which need to be 
taken into account: 

Cii ^Gn =Gn + 2Gi2lnX (C.5) 
G22 — ^ G22 = G22 + 3G23lnX (^.6) 
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which imply the following modifications at orders as 

(L) ^ (t(i) (L) + AGnL + AC, , (C.7) 

and aj: 

+ (G'nACi + AG'nCi + AGnACOL. (C.8) 

We have defined ACi = Ci — C and similarly for other quantities. This expansion is 
valid in the case where Ci multiplies qiV^Q"^) rather than q{Q'^) as was used in [p 



However the analogous expressions are straightforward to determine for that case too. 

The above formulae follow directly from the space resummed result but their 
re-interpretation, where required, as matrix projections in x space is rather straight- 
forward. We note that these formulae are currently not applicable to the improved 
broadening resummations, which have further terms that need to be taken into ac- 
count. 



D Accuracy checks 

In deriving the 'improved' resummation in section |^, the stated aim was that the result 
should be correct to within corrections of relative order a^. This requires a careful 
study of contributions that have been neglected. 

There are three main potential sources of inaccuracy that must be considered: 

• The choice of the expansion point. 

• The terms to be kept in the expansion of (and whether they need to be kept 
in the exponent). 

• The choice of terms needed in evaluating the derivative of In a, required for the 
inverse i/-transform. 



D.l Choice of expansion point 

In the standard approach of section || the radiators are both expanded around the point 
V. If one's only interest were to obtain a convergent answer for R' >2 then this could 
be achieved by keeping that expansion point and simply including the second order 
expansion of R in the exponent. 

However for R' > 2 the saddle-point of the integral is far from bQ ~ u. Indeed it is 
in a region where aslnu/bQ ~ 1, and therefore in an expansion 

R{bQ) = R{u) + ^ (\n (D.l) 

n=l ^ ^ ^ 
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where R^"'\u) is of the form ag~^ f (agL) , it is necessary to keep all terms. One way this 
could be achieved, is by not expanding in the first place and integrating, say numerically, 
over b with the full function for R{b). However this would lead to problems because for 
sufficiently large values of b one reaches a singularity associated with the Landau pole 
in R. 

If on the other hand we expand close to the saddle-point of the integral the situation 
simplifies. The width of the integrand around the saddle-point is roughly of order 
I/VR/' ~ «s"^''^) allowing us to keep a fixed number of terms in our expansion. 

For simplicity we choose to expand not around the actual saddle point, but a point 
bo close to it, differing from it by a pure numerical factor. Since in any case we will 
then have to keep terms in our expansion so as to give an accurate representation of our 
function up to ln6/6o ~ cts , the difference of a pure factor between the expansion 
and saddle points makes no difference. 



D.2 Choice of terms to be kept 

When examining the choice of terms to be kept, the discussion can be kept simpler in 
the relevant region, asL ~ 1, by noting that any quantity which is formally a^f^a^L) 
(with / some arbitrary function) is just of order a". Accordingly we will refer to R" as 
being of order ag, R'" as being of order a1 and so on. 

Keeping the first and second order terms in the expansion of TZr, our basic integral 
is of the form 

_ f dy f 1 + y/l + l/A (i„|_^)H'_l(ln|-f)2fi" .^ 

J y{i + y'r/'[ ^ J ■ ^ ^ ^ 

We see that for ~ 2 it is the R" term which ensures the integral's convergence. It 
must therefore be kept in the exponent. Examining the integral one sees that there are 
actually three possible regimes. 

• 2 — R' ^ V^s- In this case the expansion point remains close to y = 1 and the 
integral converges in a region of A In y of order 1 . 

• 2 — i?' ~ ^/<ys■ In this case the saddle and expansion points remain close toy = 1. 

— 1/2 

The relevant integration region extends down to Iny ~ —ag , but convergence 
is rapid for y > 1. 

• R' — 2 ^ ^/tts- Here —i S> y/a^ and accordingly the relevant part of the integral 
is entirely contained in the region ?/ ^ 1 and all occurrences of 1 + can simply 
be replaced by 2. 

The first region is simply that addressed in section |^, and one can neglect even the R" 
term. 
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The second and third regions require more care. We need to understand what 
happens when we multiply the integrand by a term a™(ln?/ — £)", 



^mn- J y ^2)3/2 ^ 4 j ^ «S (^iH ^ 

(D.3) 

Let us first consider the third region, where this simplifies to 

^2+R'^2i I ^^-(ln|-.)(ii'-2)-l(l„|-.)^R" . ^n. (^j^ | - ^) (D.4) 

If we are sufficiently into this region that we can neglect the {R' — 2) term then we 
obtain 

\ C m—n/2 

^mn I ttg n even 

) n odd ^ ■ ' 

The largest such contributions will come from terms such as R'"^, R^^\ and will 

all be of relative order ag, and so negligible. 

In the situation where 2 — i?' ~ the situation is more complex because the odd-n 
terms do not give zero. The reason is that the integral extends only to one side of the 
saddle-point so one loses the cancellation between lny/2 — ^ > and \ny/2 — I < 0. 
Accordingly, in this region 

- ■ , n > . 



A 

Accordingly we must keep all terms as{\ny/2 — i) and aj{lny/2 — £)^, since they con- 
tribute at the relative O (a/os) level. This is the motivation behind the set of terms 
kept in eq. ( p.6|) . 



D.3 Terms to be kept in the derivative of A 

In evaluating the inverse Mellin transform with respect to u, it is necessary to calculate 
the factor 

In the second and third regions discussed above A varies rapidly, so its derivative can 
contribute significantly to this factor and should not be neglected. 

The derivative of the full In A is technically quite complicated to evaluate because 
of the presence of the anomalous- dimension matrices. However these terms, which are 
of order as(ln?//2 — £), have two important features: they give a contribution of relative 
order ^/a^, and that contribution is relevant (and varies significantly) only in the region 
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— 1/2 ~ 

of AL ~ as . Consequently in In A when taking the derivative with respect to L, 
these extra pieces become of order ag rather than and so can be neglected. The 
same is true of all terms which contribute a relative amount y/a^ in this limited region 
of L. 

The fact that we can throw away terms contributing a relative amount ^/a^ in that 
region also allows us to modify the large-y structure of the integrand, and therefore we 
can use rather than A„. That the large-y region contributes an amount of relative 

— 1/2 

order follows from the fact that the whole integral is of order as , while the 
large-y region is unenhanced and contributes an amount of order 1. 



E Modified matching 

As was discussed in section ^, it is important that the matching respect certain proper- 
ties concerning the behaviour at the maximum of the distribution V = V^ax- We recall 
that these were: first the integrated cross section should go to whatever the correct 
upper limit happens to be, without leftover terms of order ag or higher. Secondly if 
the fixed order distribution goes to zero smoothly at the upper limit, so should the 
matched-resummed one. In e~^e~ this was always the case, whereas for many DIS vari- 
ables the distribution is non-zero at the upper limit (and the upper limit is the same 
for all orders). Here we discuss a modified matching procedure that makes sure that 
the final answer has the required behaviour in the above respects. 

The first element of the modification is to replaced 



V p 



V \Vr 



+ 1 



(E.l) 



where we have generalised through the inclusion of the power p the modification orig- 
inally proposed in (which used p = 1). Typically one might expect to consider 
1 < P ^ 3, where the upper limit is fairly arbitrary and the lower one comes from an 
assumption that the cross section contains no terms of the form l^^lnl^ with p < 1. In 
cases where we use a rescaled variable then we have 



In ^ L = - In 

XV p 



XV \ XVr 



(E.2) 



We shall also need a factor with the property that it goes to 1 rapidly for V ^ and 
to zero for V = Knax- We adopt the following form for it: 

Z{V) = l-(-^Y. (E.3) 

V •'max / 

For M and M2 matching, the replacement of L with L is usually sufficient to ensure 
that cross section goes exactly to the O (a^) upper limit. This is because at = \4iax, 

■^Certain resummations are defined, by default, not in terms of In 1 /V but rather in terms of In Vq/T^ . 
One such example is the C-parameter with Vq = 6, both in e+e^ Q and in DIS [|l3|. We write our 
formulae so that they are valid in these cases too. 
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L — 0, and accordingly ar contains no terms higher than 0{as), and the matching 
adjusts both the order ctg and terms to be exact, without introducing any additional 
terms. There is an exception to this rule for the improved broadening formula, which 
contains O {af) terms at Knax even after the replacement L ^ L — these then need to 
be subtracted. 

(2) 

If (Te (Knax) IS uou-zero (as it usually is), then the matched distribution docs not 
go to zero even if the fixed order one does, because of the alai^\Vrai,^)Qii contribution. 
Technically the most straightforward solution to the above two problems is to use the 
following formula 

a{V) = ar+[as {4'^ - ai'^) + al (^f - - alZ{V) (a^^) - 4'^) {L^G,, + L^n) 

+ 1 + asai'\V^^) + alai^\V^^) - ar{V^^)] E^(^) , (E.4) 

though strictly speaking only the Gu part of E need be raised to the power Z{V). 
Correspondingly for M2 matching we have 

a{V) = a, + as {4'^ - ai'^) + [a! (^f ' - 4'^) 

+1 + asai'\V^^) + alai^\V^^) - ar{V^^)] S^(^) . (E.5) 

For Ini? matching the situation is different in that, following the replacement L — > 
L, the matched distribution goes to zero if the fixed-order one does, but that now 
we need to fix up the value of the cross section at the upper limit. In this respect 
the procedure differs from the e~^e~ case, where replacing L ^ L led to all required 
properties automatically being satisfied. There are two reasons for the difference: the 
fact that we keep the full cr^g in front of the exponential (in e+e" the constant part is 
left out) and the fact that the fixed order cross section does not go to exactly 1 at the 
upper limit, but to 1 + C (cts). To ensure that we get exactly the same answer as the 
O («!) answer we need to insert an extra factor F in front of the exponential: 

with 

1 + CKsCTrq^ + Oi'^drf 

where have exploited the fact that o"rg''(Knax) is zero. Since F = l + O (ai?) the inclusion 
of the factor F makes no difference (at our accuracy) to the small-V behaviour of the 
answer. 

To see that the matched distribution goes to zero if the fixed-order one does, we 
restrict ourselves to the case of variables without any variable dependent scale in the 
structure function, because it turns out that none of the variables involving V in the 



• (E.7) 



36 



scale of the distribution have fixed-order distributions which go to zero at the upper 
hmit. So with this proviso one recalls 

(Trq = (1 + «^Ci)e^^i("^^)+^^("^^) Ci = . (E.8) 

Q 

In the exponent, after matching, the as and a| will have vanishing derivatives by 
construction {i.e. from the inclusion of the fixed order pieces). At order we get 
contribution only from gi and (72, but they go as and respectively and accordingly 
their derivatives vanish. This ensures that the quark resummation part gives a vanishing 
distribution. 

We also need to show that the gluon resummation piece gives a vanishing distribu- 
tion — this follows since, because we have included no g2 term in the gluon resummation 
exponent, the lowest power of L that is present is L^, automatically giving zero deriva- 
tive. 



F Parton distributions 



A complication in the practical computation of the broadening (and also Tze and r^g) 
distribution in DIS arises from the fact that the main result (|2.26|) involves operators in 
X and fiavour space, a consequence of the presence of the anomalous dimension matrix. 
In, in R^- 

The action of is of course just to change the scale of the parton distributions to 
iyQYi z.e. we can rewrite (|2.26| ) as 



T.n{B) = c 



0,7V 



V{l + 2R' 



.g-2iJ(l/B)-i?'(ln2+27B)q^^^2g2 



(F.l) 



Parton distributions are available in tabulated form as a function of scale from groups 
such as MRST CTEQ or GRV [^, so we could just choose to use these 
tabulated distributions at scale V^Q^. This is what is usually done in the context for 
example of Drell-Yan pt resummations (for recent examples, see and references 
therein). 

Instead however we have chosen to take a seemingly more complicated route and 
develop our own software for the evolution of parton distributions. There are four main 
reasons for this. 



Matching. When we carry out matching to fixed-order calculations we need to know 
quantities such as 

P(^)®g, P(^)®g, pW®pW®g, (F.2) 

where P*-^-* and P^^^ are the matrices respectively of leading and next-to-leading order 
splitting functions. By taking numerical derivatives it would be possible to obtain 
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certain combinations of the above quantities, for example 

(asP^^^ + alP^^'>^ ® q (F.3) 

from the first derivative. However this combination would be obtained only for a par- 
ticular value of ttg (that used in the original evolution). The double convolution term 
can be obtained from the second derivative of the structure functions, but only in a 
form 'polluted' by additional O (a^) terms. Furthermore some current tabulated (global 
fit) parton distribution sets use quite poor interpolation which means that numerically 
determined derivatives (especially higher derivatives) are nonsensical. 

So we in any case need to write software to calculate the quantities in ( [F.2|) . Once 
this is done, writing a full PDF evolution program involves relatively little extra work. 



Flexibility. If we follow the philosophy of e^e~ event-shape resummations, i.e. we 
include LL and NLL terms but no higher-order contributions (other than as introduced 
in the matching), then we must take eq. ( |2.26| ) literally, using only the leading-order 
splitting functions. If on the other hand we use (|F.1|) with tabulated global-fit parton 
distributions then the evolution embodied in qiV^Q^) will automatically include NLO 
splitting functions (either way we have to use NLO parton densities since we match to 
O (ctg) fixed order calculations). 

Using our own evolution code gives us a certain flexibility, and we are free to take 
a NLO parton distribution at scale Q, and then apply to it the operator exp(— _R^(_B)), 
i.e. carry out leading order evolution to scale VQ. 



Fitting CKg. When fitting for ag, for each new value of as that one wishes to examine, 
the formally correct procedure is to reevaluate all quantities using parton distributions 
fitted with that value of as. 

With the tools that are currently available for calculating fixed-order distributions, 
this means rerunning DISENT or D1SASTER++, a process which can use several tens 
of days of computing time on a modern workstation. As a result it is common to fit for 
as using a single PDF set and then to check that the results do not change significantly 
with a set corresponding to a different value of as. 

Within this approach, if we insert tabulated global-fit distributions into ( [b'.H ), then 
while most of the calculation will be done with the value of as that is explicitly inserted 
into the formulas, the single logs associated with the anomalous dimension will be 
evaluated with a value of as corresponding to the PDF set. We can study the extent 
to which this is a problem by examining how the ratio 

g(a;,Q2) 

depends on different approaches used for its evaluation (we recall our earlier definition 
for q, eq. (|4.4| )). Suppose we are calculating the broadening distribution with as(M^) = 
0.1225. The correct procedure is to calculate the ratio ( [b'.4] ) using the MRST99 hi-as 
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{as{Mz) = 0.1225) set evolved down from Q to VQ with this appropriate value of 
as{Mz)- We want to examine the error that is made if one calculates this ratio in two 
different ways (as would be done in a fit): (a) with the central-ag MRST99 (originally 
fitted with as{Mz) = 0.1175) set evaluated at scale and evolved down to V^Q^ with 
as{Mz) = 0.1175 — this is roughly equivalent to using the 'central' tabulated global-fit 
distribution in both the numerator and denominator of ([F.l|) ; (b) with the central-cts 
MRST99 set evaluated at scale and evolved down to V^Q^ with as{Mz) = 0.1225 — 
this is the philosophy of treating the anomalous dimension terms on the same footing 
as all the other single logs. The error that arises in these two different approaches is 
shown in figure ^, where one sees that the ctg used in the evolution is considerably more 
important than the used in the fit for the parton distributions. 

It should however be noted that for small-x and small-Q we can have the opposite 
situation because of the strong correlation between the value of as and the fitted gluon 
distribution. Nevertheless this analysis indicates that at the very least we want the 
option of doing our own evolution of the parton distribution. 
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Figure 5: The relative error on the ratio (|F.4|) introduced by (improper) use of the 
central MRST PDF set to calculate quantities at as{Mz) = 0.1225 (where one should 
formally use the high-a^ MRST distribution). Shown for x = 0.056, Q = 36.7 GeV 
and X = 0.18, Q = 81.3 GeV. The value of as is that used to evolve the central MRST 
parton distribution down from scale to scale V^Q^. 



Technical limitations of tabulated global-fit distributions. One final motiva- 
tion for using our evolution in calculating is qiV^Q"^) is smoothness. Generally global-fit 
parton distributions are provided in tabulated form together with an interpolating pro- 
gram. 
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In the case of the CTEQ and GRV distributions the interpolation is of reasonable 
quality. However for the current publicly available MRST distributions (dating from 
1999, [^]) only linear interpolation is used, leading to non-smoothness in Q. When 
calculating a resummed distribution (as opposed to the integrated cross section) one 
takes the derivative of ([F.l|) and this non-smoothness gets promoted to discontinuities. 
This is illustrated in figure ^ which shows the broadening distribution for Q = 81.3 GeV, 
determined in two ways: in one case we have used our own evolution between scales Q 
and VQ; in the other case we have used the MRST99 tabulated distributions to obtain 
q{V^Q'^). The clear discrepancies between the two curves (at the level of a few percent) 
are a consequence of the non-smoothness of the tabulated distribution. 

which 



We have also tested a preliminary version of the MRST 2001 distributions |39 



use improved interpolation code, and there we find that the problem is eliminated. 
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Figure 6: The broadening distribution (without power correction) calculated using 
our own (DS) evolution and with the MRST99 tabulated evolution. Shown for x = 
0.18, Q = 81.3 GeV. For other x, Q points of interest at HERA, the non-smoothness of 
the MRST distributions has a smaller effect. A leading factor of has been included 
so as to highlight the large-S region, where the problem is more severe. 



F.l Convolution and evolution algorithms. 

We have seen above that there are several motivations for evolving the parton distri- 
butions independently from the original global-fit tabulations. Accordingly we have 
developed our own evolution and convolution code. Various requirements arise from 
the need to use it for event shape resummations (though we envisage that it may well 
have wider applications): 
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• Flexibility: it should be straightforward to implement new kernels {i.e. without 
any special analytical work), since each event shape involves a new constant piece 
Ci- This also makes it straightforward to extend it to NNLL evolution. 

• Reasonable tradeoff between speed and accuracy: both setup and evolution should 
be relatively quick. 

• Robustness: we will be evolving tabulated distributions, which may have imper- 
fections due to interpolation, and which are not usually available in analytic form. 
Also, given that one may at some stage wish to use the code in a very different 
context {e.g. evolution of fragmentation functions to scales of interest for high- 
energy cosmic rays we prefer to avoid any reliance on certain properties of 
'nice behaviour' of the structure functions and their derivatives as in may re- 
duce the robustness of the algorithm (though it can enable a considerable increase 
in speed). 

The code has been written with a modular, semi object-oriented approach, within 
the hmits of the programming language used, Fortran 90. It uses an optimised x-space 
algorithm (as opposed to A^-space, where is the Mellin transform variable conjugate 
to x). Parton distributions are represented on a grid with n points uniformly spaced in 
In 1/x: 

q{x) ^ qi = q{xi) (F.5) 

where Xi = exp{—i6) with 6 the grid spacing. The parton distribution at arbitrary 
X is then defined to be equal to a linear combination of the parton distribution on 
neighbouring grid points 

i{x)+p 
j=i(x) 

where i{x) is a grid point close to x. The coefficients Cj{x) are chosen so that q{x) 
corresponds to the p^^ order interpolation of the points i{x) to i{x) +p. If one chooses 



i{xe-^) = i{x) + 1 , (F.7) 



then one has the property 



c,{x) = c,+i {xe-') . (F.8) 
With this representation of the parton distribution, convolutions then become sums: 

f{xlz)q{z) [^^q\i = 

3=0 



[P®g](a;)=/ —P{x/z)q{z)^[P®q], = Y,Pm^ (F-9) 



where 



dz 

P,,= I -P{x,/z)c,{z). (F.IO) 



z 
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Approaches of this kind have been adopted by many people, for example P2|-^, and 



indeed are a standard numerical method In the algorithms of |^,^ the grid is 
non-uniform in Inx and there are roughly 71^/2 elements of the Pij. Their evaluation 
leads to considerable overheads (analytical or numerical depending on the approach), 
while the requirement that they be held in memory during program execution (rather 
than on disk) places an upper limit on the value of n that can be used. 

As has been exploited by [HI, 331, on a uniform grid, away from the edges of the 



integral, the property (|F.8|) allows a further simplification: 

P,, = P(,+i)0+i) , (F.ll) 

i.e. Pij is only a function oii — j. Thus one has only n elements to evaluate and store, 
making it feasible to go to large values of n. Furthermore the operators for multiple 
convolutions are straightforward to evaluate, for example 

Pl = Y,P^JPJk (F.12) 

3 

and Pj| can be completely determined with O {n^) operations. In contrast if (^TTJ) 
does not hold then one needs O {n^) operations. 



In this method has been implemented only for p = 1, while has applied 
it to the general p case. For the general-p case some subtleties arise with grid edges, 
because eqs. ( [I''.6D and ( [b'.?] ) taken together imply a sum over points outside one's grid. 



The approach in |^ is to sacrifice the formal accuracy of the approach at large x, so 
that for the integration region between grid points j — 1 and j, the order is min(p, j)*'^ 
order. 

Here we observe that by choosing the i{x) appropriately (and differently according 
to the value of x used in (|F.9|) ) one can maintain the full p**^ order accuracy of the 
method. The price that we pay is that Pij now requires the computation of O (pn) 
entries rather than n, the extra entries being used to treat correctly the edge-region 
close to 2; = 1. Similarly the evaluation of an operator such as P^ requires O {pn"^) 
operations. This therefore remains much more manageable than the O (n^) factor that 
is relevant with a non-uniform grid, while preserving formal p^^ order accuracy at all 
values of x. We of course maintain the property that the evaluation of P ^ q requires 
0{n'^/2) operations 

We point out that the approach used here differs significantly in philosophy from that 
in say [^] in that the evolution in is performed with a Runge Kutta algorithm, rather 
than by a formal analytic solution to the evolution equations expressed in terms of a 
power series of P. This gives considerable simplicity because, for example, the inclusion 
of NNLL splitting functions and 3-loop running for requires no extra analytical 
calculations. 



42 



en 

Q 



< 



o) 0.95 



Q 



< 



Q 



< 




1.05 



w 0.95 

0.9 
1.1 


















^ '^'-<^-' \ ^ ^ / ^ 

_ A _ s.. V./ _ 1 > / . 


.,^.,„^ r- 1" r> /• j> r\/^ 


~i 1^ n ~i "i 'i 'i ^ "i 1 iT II 1 




MRST991 








CTEQ5M1 

CTEQ5M 

1 








(b) 





1.05 
1 



w 0.95 



0.9 



DIS scheme 




MRST99DIS1 
CTEQ5D 
CTEQ5D (MSBar evolution) 



(c) 



10" 



10" 



10" 



10 



-2 



10" 



10^ 



Figure 7: Comparisons between evolution encoded in two widely used tabulated global- 
fit PDF sets (TAB) ||2^, ^ and our evolution (DS). In all cases we start with the 
global- fit (tabulated, interpolated) distribution at 10 GeV, evolve up to 179 GeV and 
then compare to the tabulated distribution at the higher scale; (a) and (b) are in the 
MS scheme, while (c) is for the DIS scheme. 
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F.2 Testing the evolution: comparisons to tabulated global-fit 
PDF sets 



As a test of our evolution code we decided to compare to the latest available tabulated 
evolution of the MRST CTEQ |^ and GRV groups. It was our expectation 



that following the work of the 1995/96 HERA workshop (see for example [J^]), the 
latest available tabulated distributions would embody reasonably accurate evolution. 

The comparisons indicate that our code is functioning properly since we generally 
agree with the evolution of at least one of the tabulated global-fit PDF sets. However 
in certain situations (probably of limited practical phenomenological relevance) we do 
see significant disagreement with one of MRST or CTEQ and in what follows we give 
a summary of our findings. 

To carry out the comparison we adopted the following procedure. We started with 
the MRST or CTEQ tabulated (interpolated in x) distribution at Q = 10 GeV and 
evolved it up to 179 GeV with our evolution routines. We then took the ratio of the 
tabulated distribution at 179 GeV and our evolved distribution. This is plotted for 
MRST and CTEQ in figure |^. We have chosen an upper limit of 179 GeV rather than a 
'round number' because it is a grid point of the MRST distributions: the interpolation 
provided with the MRST99 PDF sets is linear in (as mentioned earlier) and a choice 
of Q in between grid points introduces errors of the order of a couple of percent. 

Figure ^ shows the results for the gluon distribution — there is good agreement 
with CTEQ and MRST over most of the range in x. However at small x < 10~^ we find 
significant disagreement with MRST at a level of up to 10%, whereas the agreement 
with CTEQ remains perfect. At large x > 0.5 there is slight disagreement with the 
CTEQ evolution, while agreement with MRST remains good up to x > 0.8. It should 
be kept in mind that in this region the gluon distribution is in any case very small, and 
poorly constrained experimentally. 

We note that the small irregularities of the curves arise because slight non-smoothness 
of the starting distributions (a consequence of the interpolation in x) gets amplified after 
convolutions with the plus-distributions of the splitting functions. 

Figure 0b shows the analogous results for the quark singlet distribution, 

S(x) = Yl + (F-13) 

j=u,d,s,... 

The problem in the MRST distributions at small x are present here too, though they 
are reduced by a factor of 2. At large x the agreement is uniformly good. We also show 
a comparison with the CTEQ5M evolution which uses some numerical approximations 
in the NLO evolution (see the note in the archive v3 of p4[)' which is useful for reference 
below. 

We have not explicitly shown comparisons with the GRV98 ms global-fit evolu- 



tion — there we find only small differences, generally less than 1%, which can 
be ascribed to a different convention for the treatment of higher-order terms in the 
evolution. Detailed comparisons with Vogt's code, using the same convention for the 
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higher-order terms, gives systematic agreement to a relative accuracy of better than 
10"^ for x < 0.8 |4l. 



We also wished to test our DIS scheme evolution, figure |^. The DIS factorisation 
scheme is defined as the scheme in which all higher-order corrections to the F2 coeffi- 
cient functions are zero. Again we start with tabulated global-fit distributions (now in 
the DIS scheme) at 10 GeV, evolve them up to 179 GeV (with DIS scheme splitting 
functions) and compare with the tabulated distributions at the higher scale. Ignoring 
the small-x problem, with MRST there is generally good agreement. There is a mod- 
erate discrepancy at large x, however we believe that this is because of the convention 
adopted by MRST to define the DIS scheme — all their fitting and evolution is done 
in the ms scheme, and then at each Q value they calculate the DIS-scheme distribution 
from the ms distributions. This differs from straightforward DIS scheme evolution by 
NNL terms, which are largest at large x because of a ln(l — a;)/(l — 2;)+ enhancement. 

With CTEQ, the agreement is quite poor and we identify two reasons for this. 
Firstly there is no updated version of the DIS scheme distributions {i.e. no equivalent 
to CTEQ5M1), so the numerical inaccuracies present in CTEQ5M (see figure |^) are 
present in CTEQ5D as well. Secondly it seems that the NLO splitting functions used 
for the evolution are in the ms scheme: if we use ms splitting functions for the evolution 
then the resulting disagreement is identical to that due to numerical approximations 
that is to be seen in CTEQ5M. 

We note that the DIS scheme was not examined in the context of the 1995/96 HERA 
workshop studies and it would perhaps be worth performing a similar 'standardis- 



ation exercise' also in the DIS scheme. 
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